### Texas in the Rearview Script 
# R version 4.3.1 (2023-06-16)

# Load Packages 
rm(list = ls()) 
library(tidyverse) 
library(ggplot2)
library(readr)
library(ggpubr)
library(sandwich)
library(lmtest)
library(stargazer)
library(plm)
library(haven)
library(dotwhisker)
library(tidyr)
library(readxl)
library(usmap)
library(dplyr)
##################################################
# Loading Primary Data 
##################################################
chairs <- read_csv("county_chairs_2019_2023.csv")
rural <- read_csv("Rural-urban Continuum Code 2013-Table 1.csv")
states <- read_csv("states.csv")
guber_2006_2022 <- read_csv("guber_2006_2022.csv")
txguber_2006_2022 <- guber_2006_2022[which(guber_2006_2022$State=="Texas"),]
texas <- read_excel("county_chairs_2008-2022.xlsx")
texas_elections <- read_csv("texaselectionclean.csv")
senate_2004_2022 <- read_csv("senate_2004_2022.csv")
uscontests <- read_excel("National Uncontested by County_Updated.xlsx")

##################################################


##################################################
# Altnernative Measures of Rural 
##################################################

rural$rural <- NA
rural$rural[rural$RUCC_2013<=3] <- 0
rural$rural[rural$RUCC_2013>=4] <- 1

rural$rurallimited <- NA
rural$rurallimited[rural$RUCC_2013<=4] <- 0
rural$rurallimited[rural$RUCC_2013>=5] <- 1

rural$rural_morelimited <- NA
rural$rural_morelimited[rural$RUCC_2013<=5] <- 0
rural$rural_morelimited[rural$RUCC_2013>=6] <- 1

rural$mostrural <- NA 
rural$mostrural[rural$RUCC_2013<=6] <- 0 
rural$mostrural[rural$RUCC_2013>=7] <- 1 
##################################################


##################################################
# Prepping Rural Identifier Data for Merge With Chair Data
##################################################

rural <- rural %>% 
  separate(County_Name, c("County_Name", "County"), " ")

rural$county <- tolower(rural$County_Name)
rural$State <- tolower(rural$State)

rural$code<- with(rural, paste0(county, State))

chairs$year <- NA 
chairs$year[chairs$`Election Cycle`=="2019-20"] <- 2020
chairs$year[chairs$`Election Cycle`=="2021-22"] <- 2022
chairs$year[chairs$`Election Cycle`=="2023-24"] <- 2024

chairs$emptychair <- NA
chairs$emptychair[chairs$`Chair Name`=="VACANT"] <- 1
chairs$emptychair[chairs$`Chair Name`!="VACANT"] <- 0


chairs <- merge(chairs, states, by.x="State", by.y = "State_Name")

chairs$county <- tolower(chairs$`County Name`)
chairs$Abbr <- tolower(chairs$Abbr)

chairs$code<- with(chairs, paste0(county, Abbr))

# Nationl Chairs Data 
dat <- merge(chairs, rural, by="code")
dat$spot <- 1


##################################################
##### Plot of Rural-Urban Democratic/Republican Vacancies Nationally 
##################################################

# Metro-non metro 
results <- dat %>%
  group_by(rural, year, Party) %>%
  summarize(vacancies = (sum(emptychair,  na.rm=T)))


results2 <- dat %>%
  group_by(rural, year, Party) %>%
  summarize(vacancies = (sum(emptychair,  na.rm=T)/sum(spot,  na.rm=T)))


results2$group <- NA 
results2$group[results2$rural==1 & results2$Party=="Democratic"] <- "Rural-Democratic"
results2$group[results2$rural==0 & results2$Party=="Democratic"] <- "Urban-Democratic"
results2$group[results2$rural==1 & results2$Party=="Republican"] <- "Rural-Republican"
results2$group[results2$rural==0 & results2$Party=="Republican"] <- "Urban-Republican"

results$group <- NA 
results$group[results$rural==1 & results$Party=="Democratic"] <- "Rural-Democratic"
results$group[results$rural==0 & results$Party=="Democratic"] <- "Urban-Democratic"
results$group[results$rural==1 & results$Party=="Republican"] <- "Rural-Republican"
results$group[results$rural==0 & results$Party=="Republican"] <- "Urban-Republican"



p1 <- ggplot(results, aes(x=year, y= vacancies, color=group, fill=group)) +   
  geom_smooth(method="loess", se=F) + 
  labs(y="Number of County Chair Vacancies", 
       x = "Year", 
       color="Geography-Party", 
       fill="Geography-Party") + scale_fill_grey() + scale_color_grey() + 
  theme_classic()
p1 

p2 <- ggplot(results2, aes(x=year, y= vacancies, color=group, fill=group)) +   
  geom_smooth(method="loess", se=F) + 
  labs(y="Proportion of County Chair Vacancies", 
       x = "Year", 
       color="Geography-Party", 
       fill="Geography-Party") + scale_fill_grey() + scale_color_grey() + 
  theme_classic()
p2 


##################################################
# Variation on "Rural` Coding for Robustness (Appendix)
##################################################

results2 <- dat %>%
  group_by(rurallimited, year, Party) %>%
  summarize(vacancies = (sum(emptychair,  na.rm=T)/sum(spot,  na.rm=T)))

results2$group <- NA 
results2$group[results2$rurallimited==1 & results2$Party=="Democratic"] <- "Rural-Democratic"
results2$group[results2$rurallimited==0 & results2$Party=="Democratic"] <- "Urban-Democratic"
results2$group[results2$rurallimited==1 & results2$Party=="Republican"] <- "Rural-Republican"
results2$group[results2$rurallimited==0 & results2$Party=="Republican"] <- "Urban-Republican"

results2 <- dat %>%
  group_by(rural_morelimited, year, Party) %>%
  summarize(vacancies = (sum(emptychair,  na.rm=T)/sum(spot,  na.rm=T)))


results2$group <- NA 
results2$group[results2$rural_morelimited==1 & results2$Party=="Democratic"] <- "Rural-Democratic"
results2$group[results2$rural_morelimited==0 & results2$Party=="Democratic"] <- "Urban-Democratic"
results2$group[results2$rural_morelimited==1 & results2$Party=="Republican"] <- "Rural-Republican"
results2$group[results2$rural_morelimited==0 & results2$Party=="Republican"] <- "Urban-Republican"

rural$rural_morelimited <- NA
rural$rural_morelimited[rural$RUCC_2013<=5] <- 0
rural$rural_morelimited[rural$RUCC_2013>=6] <- 1



results2 <- dat %>%
  group_by(mostrural, year, Party) %>%
  summarize(vacancies = (sum(emptychair,  na.rm=T)/sum(spot,  na.rm=T)))


results2$group <- NA 
results2$group[results2$mostrural==1 & results2$Party=="Democratic"] <- "Rural-Democratic"
results2$group[results2$mostrural==0 & results2$Party=="Democratic"] <- "Urban-Democratic"
results2$group[results2$mostrural==1 & results2$Party=="Republican"] <- "Rural-Republican"
results2$group[results2$mostrural==0 & results2$Party=="Republican"] <- "Urban-Republican"


dat$ruralsub <- NA 
dat$ruralsub[dat$RUCC_2013<=2] <- 0
dat$ruralsub[dat$RUCC_2013>=3] <- 1


results2 <- dat %>%
  group_by(ruralsub, year, Party) %>%
  summarize(vacancies = (sum(emptychair,  na.rm=T)/sum(spot,  na.rm=T)))


results2$group <- NA 
results2$group[results2$ruralsub==1 & results2$Party=="Democratic"] <- "Rural-Democratic"
results2$group[results2$ruralsub==0 & results2$Party=="Democratic"] <- "Urban-Democratic"
results2$group[results2$ruralsub==1 & results2$Party=="Republican"] <- "Rural-Republican"
results2$group[results2$ruralsub==0 & results2$Party=="Republican"] <- "Urban-Republican"

dat$ruralsubsub <- NA 
dat$ruralsubsub[dat$RUCC_2013<=1] <- 0
dat$ruralsubsub[dat$RUCC_2013>=2] <- 1


results2 <- dat %>%
  group_by(ruralsubsub, year, Party) %>%
  summarize(vacancies = (sum(emptychair,  na.rm=T)/sum(spot,  na.rm=T)))


results2$group <- NA 
results2$group[results2$ruralsubsub==1 & results2$Party=="Democratic"] <- "Rural-Democratic"
results2$group[results2$ruralsubsub==0 & results2$Party=="Democratic"] <- "Urban-Democratic"
results2$group[results2$ruralsubsub==1 & results2$Party=="Republican"] <- "Rural-Republican"
results2$group[results2$ruralsubsub==0 & results2$Party=="Republican"] <- "Urban-Republican"



##################################################
### Prepping Texas Chair Data 
##################################################

texas$state <- "tx"
texas$COUNTY <- tolower(texas$COUNTY)
texas$code <-with(texas, paste0(COUNTY , state))

texas <- merge(texas, rural, by="code")
texas$spot <- 1

texas$emptychair <- NA
texas$emptychair[texas$NAME=="VACANT"] <- 1
texas$emptychair[texas$NAME!="VACANT"] <- 0

texas$emptychair[is.na(texas$emptychair)] <- 0

texas_2022 <- texas[which(texas$ELECTION_YEAR==2022),]
texas_2022d <- texas_2022[which(texas_2022$PARTY=="Dems"),] #  0.2405 
texas_2022dr <- texas_2022d[which(texas_2022d$rural==1),] #  0.3082
texas_2022du <- texas_2022d[which(texas_2022d$rural==0),] #  0.3082

texas2 <- texas[which(texas$PARTY=="Dems"),]
texas2 <- texas2 %>%
  group_by(FIPS) %>%
  summarize(vacancies = (sum(emptychair,  na.rm=T)))

texas2$vacancies[texas2$vacancies>=1] <- 1

##################################################
# Making County-Level Maps and Plots for Main Text
##################################################

countypop <- usmap::countypop
countypop$FIPS <- as.numeric(countypop$fips)

countypop2 <- merge(countypop, texas2, by="FIPS")

countypop2[is.na(countypop2$emptychair),] <- 0

# Texas Map Plot Main Text
plot_usmap(data = countypop2 , values = "vacancies", include = c("TX")) + 
  scale_fill_continuous(
    low = "white", high = "black", na.value="white", label = scales::comma)  + labs(fill="Vacant Democratic Chair Ever")  + theme(legend.position="none")


# Texas Vacancies Over Time Plot Main Text
results <- texas%>%
  group_by(rural, ELECTION_YEAR, PARTY) %>%
  summarize(vacancies = (sum(emptychair,  na.rm=T)/sum(spot,  na.rm=T)))

results$group <- NA 
results$group[results$rural==1 & results$PARTY=="Dems"] <- "Rural-Democratic"
results$group[results$rural==0 & results$PARTY=="Dems"] <- "Urban-Democratic"
results$group[results$rural==1 & results$PARTY=="Reps"] <- "Rural-Republican"
results$group[results$rural==0 & results$PARTY=="Reps"] <- "Urban-Republican"

ggplot(results, aes(x=ELECTION_YEAR, y= vacancies, color=group, fill=group)) +   
  geom_smooth(method="loess", se=F) + 
  labs(y="Number of County Chair Vacancies (Texas)", 
       x = "Year", 
       color="Geography-Party", 
       fill="Geography-Party") + scale_fill_grey() + scale_color_grey() + 
  theme_classic()


# Alternative Measures of Rural Vacancies
results2 <- texas%>%
  group_by(rural_morelimited, ELECTION_YEAR, PARTY) %>%
  summarize(vacancies = (sum(emptychair,  na.rm=T)/sum(spot,  na.rm=T)))

texas$mostrural <- NA 
texas$mostrural[texas$RUCC_2013<=6] <- 0 
texas$mostrural[texas$RUCC_2013>=7] <- 1 

results2 <- texas%>%
  group_by(mostrural, ELECTION_YEAR, PARTY) %>%
  summarize(vacancies = (sum(emptychair,  na.rm=T)/sum(spot,  na.rm=T)))


texas$ruralsub <- NA 
texas$ruralsub[texas$RUCC_2013<=3] <- 0 
texas$ruralsub[texas$RUCC_2013>=4] <- 1 

results2 <- texas%>%
  group_by(ruralsub, ELECTION_YEAR, PARTY) %>%
  summarize(vacancies = (sum(emptychair,  na.rm=T)/sum(spot,  na.rm=T)))


##################################################
# County Level Elections Analyses 
##################################################

# Texas County Presidential Vote Share Analysis 
toy1 <- lm(demshare ~ emptydem + as.factor(FIPS) + as.factor(ELECTION_YEAR), data=texas_elections, w=totalvotes)
toy2 <- lm(demshare ~ emptygop + as.factor(FIPS) + as.factor(ELECTION_YEAR), data=texas_elections, w=totalvotes)

toy1cl <- coeftest(toy1, vcov = vcovCL, cluster = ~FIPS)
toy2cl <- coeftest(toy2, vcov = vcovCL, cluster = ~FIPS)


stargazer(toy1, toy2,
          dep.var.labels=c("Democratic Two-Party Vote Share"),
          covariate.labels=c("No Democratic Chair", "No Republican Chair") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 

stargazer(toy1cl, toy2cl,
          dep.var.labels=c("Democratic Two-Party Vote Share"),
          covariate.labels=c("No Democratic Chair", "No Republican Chair") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 


## Texas County Election Robustness with Lagged Dem Vote Share & Rural FEs (Appendix)
texas_elections2 <- 
  texas_elections %>%
  group_by(FIPS) %>%
  mutate(demsharelag = dplyr::lag(demshare, n = 1, default = NA))

# Results Controlling For Previous Dem Share (and allowing for rural trend)
toy1 <- lm(demshare ~ emptydem + demsharelag + as.factor(FIPS) + as.factor(ELECTION_YEAR), data=texas_elections2, w=totalvotes)
toy2 <- lm(demshare ~ emptygop + demsharelag + as.factor(FIPS) + as.factor(ELECTION_YEAR), data=texas_elections2, w=totalvotes)
toy3 <- lm(demshare ~ emptydem + demsharelag + rural + as.factor(FIPS) + as.factor(ELECTION_YEAR), data=texas_elections2, w=totalvotes)

toy1cl <- coeftest(toy1, vcov = vcovCL, cluster = ~FIPS)
toy2cl <- coeftest(toy2, vcov = vcovCL, cluster = ~FIPS)
toy3cl <- coeftest(toy3, vcov = vcovCL, cluster = ~FIPS)

# Panel with lagged DV 
stargazer(toy1, toy2, toy3, 
          dep.var.labels=c("Democratic Two-Party Vote Share"),
          covariate.labels=c("No Democratic Chair", "No Republican Chair") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 


# Panel with lagged DV 
stargazer(toy1cl, toy2cl,toy3cl,
          dep.var.labels=c("Democratic Two-Party Vote Share"),
          covariate.labels=c("No Democratic Chair", "No Republican Chair") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 


texas_elections2 <- merge(texas_elections2, rural, by="FIPS")


# Turnout Is Estimated based on total population 
texas_elections2$turnout_est <- texas_elections2$totalvotes/texas_elections2$Population_2010

## TEXAS COUNTY-LEVEL TURNOUT ANALYSIS
toy1 <- lm(turnout_est ~ emptydem + as.factor(FIPS) + as.factor(ELECTION_YEAR), data=texas_elections2, w=totalvotes)
toy2 <- lm(turnout_est ~ emptygop  + as.factor(FIPS) + as.factor(ELECTION_YEAR), data=texas_elections2, w=totalvotes)

stargazer(toy1, toy2, 
          dep.var.labels=c("Voter Turnout"),
          covariate.labels=c("No Democratic Chair", "No Republican Chair") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 

toy1cl <- coeftest(toy1, vcov = vcovCL, cluster = ~FIPS)
toy2cl <- coeftest(toy2, vcov = vcovCL, cluster = ~FIPS)

stargazer(toy1cl, toy2cl,
          dep.var.labels=c("Voter Turnout"),
          covariate.labels=c("No Democratic Chair", "No Republican Chair") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 


## TEXAS COUNTY-LEVEL TURNOUT ANALYSIS (WITH RURAL FE)
toy1 <- lm(turnout_est ~ emptydem  + rural.x + as.factor(FIPS) + as.factor(ELECTION_YEAR), data=texas_elections2, w=totalvotes)
toy2 <- lm(turnout_est ~ emptygop + rural.x + as.factor(FIPS) + as.factor(ELECTION_YEAR), data=texas_elections2, w=totalvotes)

stargazer(toy1, toy2, 
          dep.var.labels=c("Voter Turnout"),
          covariate.labels=c("No Democratic Chair", "No Republican Chair") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 

toy1cl <- coeftest(toy1, vcov = vcovCL, cluster = ~FIPS)
toy2cl <- coeftest(toy2, vcov = vcovCL, cluster = ~FIPS)

stargazer(toy1cl, toy2cl,
          dep.var.labels=c("Voter Turnout"),
          covariate.labels=c("No Democratic Chair", "No Republican Chair") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 


########## Senate Elections County-Level Results
senate_2004_2022 <- merge(senate_2004_2022, states, by.x = "State", by.y = "State_Name")
senate_2004_2022$Area <- tolower(senate_2004_2022$Area)
senate_2004_2022$Abbr <- tolower(senate_2004_2022$Abbr)

senate_2004_2022$code <- with(senate_2004_2022, paste0(Area, Abbr))
senate_2004_2022 <- merge(senate_2004_2022, rural, by="code")
senate_2004_2022_tx <- merge(senate_2004_2022, texas_wide, by=c("FIPS", "ELECTION_YEAR"))# merge like this for nationwide
senate_2004_2022_tx$DemVotesMajorPercent <- senate_2004_2022_tx$DemVotesMajorPercent/100

toy1 <- lm(DemVotesMajorPercent ~ emptydem + as.factor(FIPS) + as.factor(ELECTION_YEAR), data=senate_2004_2022_tx, w=Population_2010)
toy2 <- lm(DemVotesMajorPercent ~ emptydem + rural.x + as.factor(FIPS) + as.factor(ELECTION_YEAR), data=senate_2004_2022_tx, w=Population_2010)


stargazer(toy1, toy2, 
          dep.var.labels=c("Democratic Senate Vote"),
          covariate.labels=c("No Democratic Chair", "Rural") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 

toy1cl <- coeftest(toy1, vcov = vcovCL, cluster = ~FIPS)
toy2cl <- coeftest(toy2, vcov = vcovCL, cluster = ~FIPS)

stargazer(toy1cl, toy2cl,
          dep.var.labels=c("Democratic Senate Vote"),
          covariate.labels=c("No Democratic Chair", "Rural") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 



######  Gubernatorial Elections County-Level Results
texas_widet$Year <- texas_widet$ELECTION_YEAR
texas_widet$Area <- toupper(texas_widet$County_Name)
txguber <- merge(texas_widet, txguber_2006_2022, by=c("Area", "Year"))

txguber$DemVotesMajorPercent <- txguber$DemVotesMajorPercent/100
toy1 <- lm(DemVotesMajorPercent ~ emptydem + as.factor(FIPS) + as.factor(ELECTION_YEAR), data=txguber, w=Population_2010)
toy2 <- lm(DemVotesMajorPercent ~ emptydem + rural.x + as.factor(FIPS) + as.factor(ELECTION_YEAR), data=txguber, w=Population_2010)


stargazer(toy1, toy2, 
          dep.var.labels=c("Democratic Gubernatorial Vote"),
          covariate.labels=c("No Democratic Chair", "Rural") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 

toy1cl <- coeftest(toy1, vcov = vcovCL, cluster = ~FIPS)
toy2cl <- coeftest(toy2, vcov = vcovCL, cluster = ~FIPS)

stargazer(toy1cl, toy2cl,
          dep.var.labels=c("Democratic Gubernatorial Vote"),
          covariate.labels=c("No Democratic Chair", "Rural") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 





##################################################
#   CES Analyses of Chair Vacancies
##################################################
ces <- read_dta("cumulative_2006-2023.dta")
CES_CLEANED <- read_csv("CES_CLEANED.csv")
CESTEXAS <- read_csv("cestexas.csv")

CES_CLEANED <- merge(ces, CES_CLEANED,  by.y="caseid", by.x="case_id", all=T)
CES_CLEANED2 <- merge(CES_CLEANED, rural, by="FIPS")



##################################################
# Plotting National Chair Vacancy data 
##################################################
keep <- c("year", "FIPS", "rural", "emptychair", "Party", "code", "County_Name")

dat2 <- dat[keep]

dat2_wide <- spread(dat2, Party,emptychair)

dat2_wide$emptydem <- dat2_wide$Democratic
dat2_wide$emptygop <- dat2_wide$Republican

dat2_wide$spot <- 1
dat2_wide <- merge(dat2_wide, rural, by="FIPS")

dat2_wide2 <- dat2_wide[which(dat2_wide$year==2022),]



# Vacancies per state in 2022 
states  <- dat2_wide2 %>%
  group_by(State) %>%
  summarize(vacancies = (sum(emptydem,  na.rm=T)/sum(spot,  na.rm=T)))
# Vacancies ever 
treatedever <- dat2_wide %>%
  group_by(FIPS) %>%
  summarize(treatedever = sum(emptydem,  na.rm=T))
# By Year 
treatedever_year <- dat2_wide %>%
  group_by(FIPS, year) %>%
  summarize(treatedever = sum(emptydem,  na.rm=T))

# All counties that lost a Dem chair and when
treatedever_year <- treatedever_year[which(treatedever_year$treatedever==1),]

statepop <- statepop
states$abbr <- toupper(states$State)
statepop2 <- merge(statepop, states, by="abbr", all=T)

statepop2[is.na(statepop2$vacancies),] <- 0 # coding missing data at 0 for visualization

# Map of Vacanies in US OVERALL BY STATE
plot_usmap(data = statepop2 , values = "vacancies") + 
  scale_fill_continuous(
    low = "white", high = "black", na.value="white", label = scales::comma)  + labs(fill="Vacant Democratic Chair Proportion")  + theme(legend.position="bottom")

keep <- c("year", "FIPS", "rural.x", "emptydem", "emptygop", "code.x")

dat2_wide <- dat2_wide[keep]

CES<- merge(CES_CLEANED, dat2_wide, by.x=c("FIPS", "Year"), by.y=c("FIPS", "year"))


##################################################
# Recoding CES data 
##################################################

CES$black <- 0
CES$black[CES$race==2] <- 1

CES$latino <- 0
CES$latino[CES$race==3] <- 1

CES$faminc[CES$faminc==13] <- NA

CES$republican <- NA 
CES$republican[CES$pid7==1] <- -3
CES$republican[CES$pid7==2] <- -2
CES$republican[CES$pid7==3] <- -1
CES$republican[CES$pid7==4] <- 0
CES$republican[CES$pid7==5] <-1
CES$republican[CES$pid7==6] <-2
CES$republican[CES$pid7==7] <-3

CES$female <- 0
CES$female[CES$gender==2] <-1 


CES$DEMVOTE <- NA
CES$DEMVOTE[CES$presvote20post==1] <- 1
CES$DEMVOTE[CES$presvote20post==2] <- 0

# Creating Rural Only Data Set
CESrural <- CES[which(CES$rural.x==1),]

# Coding DVS
CESrural$turnout <- NA 
CESrural$turnout[CESrural$vv_turnout_gvm==1] <- 1 
CESrural$turnout[CESrural$vv_turnout_gvm==2] <- 0 
CESrural$turnout[CESrural$vv_turnout_gvm==3] <- 0 

CESrural$voteddemsentp <- NA 
CESrural$voteddemsentp[CESrural$voted_sen_party==1] <- 1 
CESrural$voteddemsentp[CESrural$voted_sen_party==2] <- 0


CESrural$voteddemhousetp<- NA 
CESrural$voteddemhousetp[CESrural$voted_rep_party==1] <- 1 
CESrural$voteddemhousetp[CESrural$voted_rep_party==2] <- 0


CESrural$DEMVOTE <- NA
CESrural$DEMVOTE[CESrural$voted_pres_party==1] <- 1
CESrural$DEMVOTE[CESrural$voted_pres_party==2] <- 0


CESrural$votedgovdem <- NA
CESrural$votedgovdem[CESrural$voted_gov_party==1] <- 1
CESrural$votedgovdem[CESrural$voted_gov_party==2] <- 0



## CES Analyses
m1 <- lm(contacted ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr + as.factor(Year) + as.factor(state), data=CESrural, w=commonweight)
m2 <- lm(turnout ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr + as.factor(Year)+ as.factor(state), data=CESrural, w=commonweight)
m3 <- lm(DEMVOTE ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr+ as.factor(state), data=CESrural, w=commonweight)
m4 <- lm(voteddemhousetp~ emptydem + republican + educ + faminc + female + black +  latino + birthyr+ as.factor(Year)+ as.factor(state), data=CESrural, w=commonweight)
m5 <- lm(voteddemsentp ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr + as.factor(Year)+ as.factor(state), data=CESrural, w=commonweight)
m6 <- lm(votedgovdem ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr + as.factor(Year)+ as.factor(state), data=CESrural, w=commonweight)

stargazer(m1, m2, m3, m4, m5, m6, 
          dep.var.labels=c("Contacted", "Validated Turnout",  "Voted Dem Pres" ,"Voted Dem House" , "Voted Dem Senate", "Voted Dem Governor"),
          covariate.labels=c("No Democratic Chair"), style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 

toycl <- coeftest(m1, vcov = vcovCL, cluster = ~FIPS)
toy1cl <- coeftest(m2, vcov = vcovCL, cluster = ~FIPS)
toy2cl <- coeftest(m3, vcov = vcovCL, cluster = ~FIPS)
toy3cl <- coeftest(m4, vcov = vcovCL, cluster = ~FIPS)
toy4cl <- coeftest(m5, vcov = vcovCL, cluster = ~FIPS)
toy5cl <- coeftest(m6, vcov = vcovCL, cluster = ~FIPS)


stargazer(toycl, toy1cl, toy2cl, toy3cl, toy4cl, toy5cl,
          dep.var.labels=c("Contacted", "Voted Dem President" , "Voted Dem House" , "Voted Dem Senate", "Voted Dem Governor"),
          covariate.labels=c("No Democratic Chair") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 


## CES Analyses with Alternative Measures of Rural 

CESrural2 <- merge(CESrural, rural, by="FIPS")

CESrural_l1 <- CESrural2[which(CESrural2$rurallimited==1),]
m1 <- lm(contacted ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr + as.factor(Year) + as.factor(state), data=CESrural_l1, w=commonweight)
m2 <- lm(turnout ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr + as.factor(Year)+ as.factor(state), data=CESrural_l1, w=commonweight)
m3 <- lm(DEMVOTE ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr+ as.factor(state), data=CESrural_l1, w=commonweight)
m4 <- lm(voteddemhousetp~ emptydem + republican + educ + faminc + female + black +  latino + birthyr+ as.factor(Year)+ as.factor(state), data=CESrural_l1, w=commonweight)
m5 <- lm(voteddemsentp ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr + as.factor(Year)+ as.factor(state), data=CESrural_l1, w=commonweight)
m6 <- lm(votedgovdem ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr + as.factor(Year)+ as.factor(state), data=CESrural_l1, w=commonweight)


stargazer(m1, m2, m3, m4, m5, m6, 
          dep.var.labels=c("Contacted", "Validated Turnout",  "Voted Dem Pres" ,"Voted Dem House" , "Voted Dem Senate", "Voted Dem Governor"),
          covariate.labels=c("No Democratic Chair"), style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 

toycl <- coeftest(m1, vcov = vcovCL, cluster = ~FIPS)
toy1cl <- coeftest(m2, vcov = vcovCL, cluster = ~FIPS)
toy2cl <- coeftest(m3, vcov = vcovCL, cluster = ~FIPS)
toy3cl <- coeftest(m4, vcov = vcovCL, cluster = ~FIPS)
toy4cl <- coeftest(m5, vcov = vcovCL, cluster = ~FIPS)
toy5cl <- coeftest(m6, vcov = vcovCL, cluster = ~FIPS)


stargazer(toycl, toy1cl, toy2cl, toy3cl, toy4cl, toy5cl,
          dep.var.labels=c("Contacted", "Validated Turnout", "Voted Dem President" , "Voted Dem House" , "Voted Dem Senate", "Voted Dem Governor"),
          covariate.labels=c("No Democratic Chair") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 


CESrural_l2 <- CESrural2[which(CESrural2$rural_morelimited==1),]

m1 <- lm(contacted ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr + as.factor(Year) + as.factor(state), data=CESrural_l2, w=commonweight)
m2 <- lm(turnout ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr + as.factor(Year)+ as.factor(state), data=CESrural_l2, w=commonweight)
m3 <- lm(DEMVOTE ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr+ as.factor(state), data=CESrural_l2, w=commonweight)
m4 <- lm(voteddemhousetp~ emptydem + republican + educ + faminc + female + black +  latino + birthyr+ as.factor(Year)+ as.factor(state), data=CESrural_l2, w=commonweight)
m5 <- lm(voteddemsentp ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr + as.factor(Year)+ as.factor(state), data=CESrural_l2, w=commonweight)
m6 <- lm(votedgovdem ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr + as.factor(Year)+ as.factor(state), data=CESrural_l2, w=commonweight)


stargazer(m1, m2, m3, m4, m5, m6, 
          dep.var.labels=c("Contacted", "Validated Turnout",  "Voted Dem Pres" ,"Voted Dem House" , "Voted Dem Senate", "Voted Dem Governor"),
          covariate.labels=c("No Democratic Chair"), style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 

toycl <- coeftest(m1, vcov = vcovCL, cluster = ~FIPS)
toy1cl <- coeftest(m2, vcov = vcovCL, cluster = ~FIPS)
toy2cl <- coeftest(m3, vcov = vcovCL, cluster = ~FIPS)
toy3cl <- coeftest(m4, vcov = vcovCL, cluster = ~FIPS)
toy4cl <- coeftest(m5, vcov = vcovCL, cluster = ~FIPS)
toy5cl <- coeftest(m6, vcov = vcovCL, cluster = ~FIPS)


stargazer(toycl, toy1cl, toy2cl, toy3cl, toy4cl, toy5cl,
          dep.var.labels=c("Contacted", "Validated Turnout", "Voted Dem President" , "Voted Dem House" , "Voted Dem Senate", "Voted Dem Governor"),
          covariate.labels=c("No Democratic Chair") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 


CESrural_l3<- CESrural2[which(CESrural2$mostrural==1),]


m1 <- lm(contacted ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr + as.factor(Year) + as.factor(state), data=CESrural_l3, w=commonweight)
m2 <- lm(turnout ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr + as.factor(Year)+ as.factor(state), data=CESrural_l3, w=commonweight)
m3 <- lm(DEMVOTE ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr+ as.factor(state), data=CESrural_l3, w=commonweight)
m4 <- lm(voteddemhousetp~ emptydem + republican + educ + faminc + female + black +  latino + birthyr+ as.factor(Year)+ as.factor(state), data=CESrural_l3, w=commonweight)
m5 <- lm(voteddemsentp ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr + as.factor(Year)+ as.factor(state), data=CESrural_l3, w=commonweight)
m6 <- lm(votedgovdem ~ emptydem + republican + educ + faminc + female + black +  latino + birthyr + as.factor(Year)+ as.factor(state), data=CESrural_l3, w=commonweight)


stargazer(m1, m2, m3, m4, m5, m6, 
          dep.var.labels=c("Contacted", "Validated Turnout",  "Voted Dem Pres" ,"Voted Dem House" , "Voted Dem Senate", "Voted Dem Governor"),
          covariate.labels=c("No Democratic Chair"), style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 

toycl <- coeftest(m1, vcov = vcovCL, cluster = ~FIPS)
toy1cl <- coeftest(m2, vcov = vcovCL, cluster = ~FIPS)
toy2cl <- coeftest(m3, vcov = vcovCL, cluster = ~FIPS)
toy3cl <- coeftest(m4, vcov = vcovCL, cluster = ~FIPS)
toy4cl <- coeftest(m5, vcov = vcovCL, cluster = ~FIPS)
toy5cl <- coeftest(m6, vcov = vcovCL, cluster = ~FIPS)


stargazer(toycl, toy1cl, toy2cl, toy3cl, toy4cl, toy5cl,
          dep.var.labels=c("Contacted", "Validated Turnout", "Voted Dem President" , "Voted Dem House" , "Voted Dem Senate", "Voted Dem Governor"),
          covariate.labels=c("No Democratic Chair") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 




##################################################
#  National Uncontested Seats Results 
##################################################

# Prepping Data
statepop <- statepop

states <- merge(states, statepop, by="abbr")

uscontests$County <- tolower(uscontests$County)
states$State <- states$full

uscontests <- merge(uscontests, states, by="State")

uscontests$State2 <- tolower(uscontests$abbr)

uscontests$code.x <-with(uscontests, paste0(County , State2))

uscontests <- merge(uscontests, rural, by.x="code.x", by.y = "code" )

uscontests$uncontesteddem <- NA 
uscontests$uncontesteddem[uscontests$`Democrat (Y/N)`=="Y"] <- 0 
uscontests$uncontesteddem[uscontests$`Democrat (Y/N)`=="N"] <- 1

uscontests$uncontestedgop<- NA 
uscontests$uncontestedgop[uscontests$`Republican (Y/N)`=="Y"] <- 0 
uscontests$uncontestedgop[uscontests$`Republican (Y/N)`=="N"] <- 1

uscontests$spot <- 1


##################################################
# Uncontested Seats Plot
##################################################


uscontests2 <- uscontests


results <- uscontests2 %>%
  group_by(rural, Year) %>%
  summarize(uncontested  = (sum(uncontesteddem,  na.rm=T)/sum(spot,  na.rm=T)))
results$Party <- "Democratic"

results2 <- uscontests2 %>%
  group_by(rural, Year) %>%
  summarize(uncontested  = (sum(uncontestedgop,  na.rm=T)/sum(spot,  na.rm=T)))


results2$Party <- "Republican"


resultsfull <- rbind(results, results2)

resultsfull$group <- NA 
resultsfull$group[resultsfull$rural==1 & resultsfull$Party=="Democratic"] <- "Rural-Democratic"
resultsfull$group[resultsfull$rural==0 & resultsfull$Party=="Democratic"] <- "Urban-Democratic"
resultsfull$group[resultsfull$rural==1 & resultsfull$Party=="Republican"] <- "Rural-Republican"
resultsfull$group[resultsfull$rural==0 & resultsfull$Party=="Republican"] <- "Urban-Republican"


ggplot(resultsfull, aes(x=Year, y= uncontested, color=group, fill=group)) +   
  geom_line(size=0.75) + 
  labs(y="Proportion with Uncontested Congressional Seat", 
       x = "Year", 
       color="Geography-Party", 
       fill="Geography-Party") + scale_fill_grey() + scale_color_grey() + 
  theme_classic()




##################################################
#  National Uncontested Seats CES Analysis
##################################################

uscontests <- read_excel("National Uncontested by County_Updated.xlsx")

uscontests$County <- tolower(uscontests$County)

uscontests <- merge(uscontests, states, by = "State")

uscontests$State2 <- tolower(uscontests$abbr)

uscontests$code.x <-with(uscontests, paste0(County , State2))

uscontests <- merge(uscontests, rural, by.x="code.x", by.y = "code" )



uscontests$uncontesteddem <- NA 
uscontests$uncontesteddem[uscontests$`Democrat (Y/N)`=="Y"] <- 0 
uscontests$uncontesteddem[uscontests$`Democrat (Y/N)`=="N"] <- 1

uscontests$uncontestedgop<- NA 
uscontests$uncontestedgop[uscontests$`Republican (Y/N)`=="Y"] <- 0 
uscontests$uncontestedgop[uscontests$`Republican (Y/N)`=="N"] <- 1

uscontests$spot <- 1

uscontests2 <- uscontests


# Plot of Uncontested Seats Over Time 
results <- uscontests2 %>%
  group_by(rural, Year) %>%
  summarize(uncontested  = (sum(uncontesteddem,  na.rm=T)/sum(spot,  na.rm=T)))
results$Party <- "Democratic"

results2 <- uscontests2 %>%
  group_by(rural, Year) %>%
  summarize(uncontested  = (sum(uncontestedgop,  na.rm=T)/sum(spot,  na.rm=T)))

results2$Party <- "Republican"


resultsfull <- rbind(results, results2)

resultsfull$group <- NA 
resultsfull$group[resultsfull$rural==1 & resultsfull$Party=="Democratic"] <- "Rural-Democratic"
resultsfull$group[resultsfull$rural==0 & resultsfull$Party=="Democratic"] <- "Urban-Democratic"
resultsfull$group[resultsfull$rural==1 & resultsfull$Party=="Republican"] <- "Rural-Republican"
resultsfull$group[resultsfull$rural==0 & resultsfull$Party=="Republican"] <- "Urban-Republican"


ggplot(resultsfull, aes(x=Year, y= uncontested, color=group, fill=group)) +   
  geom_smooth(method="loess", se=F) + 
  labs(y="Proportion of Counties with Uncontested Congressional Seat", 
       x = "Year", 
       color="Geography-Party", 
       fill="Geography-Party") + scale_fill_grey() + scale_color_grey() + 
  theme_classic()


ggplot(resultsfull, aes(x=Year, y= uncontested, color=group, fill=group)) +   
  geom_line(size=0.75) + 
  labs(y="Proportion with Uncontested Congressional Seat", 
       x = "Year", 
       color="Geography-Party", 
       fill="Geography-Party") + scale_fill_grey() + scale_color_grey() + 
  theme_classic()


# Uncontested Seats Data Merged with CES
surveyuncon <- merge(uscontests2, CES_CLEANED2, by.x = c("Year", "FIPS"), by.y = c("Year", "FIPS"))

# Coding Variables for Analysis 
surveyuncon$female[surveyuncon$gender==1]<- 0
surveyuncon$female[surveyuncon$gender==2]<- 1

surveyuncon$black <- 0
surveyuncon$black[surveyuncon$race==2] <- 1

surveyuncon$latinx <- 0
surveyuncon$latinx[surveyuncon$race==3] <- 1

surveyuncon$republican <- NA
surveyuncon$republican[surveyuncon$pid7==5] <- 1 
surveyuncon$republican[surveyuncon$pid7==6] <- 2 
surveyuncon$republican[surveyuncon$pid7==7] <- 3
surveyuncon$republican[surveyuncon$pid7==4] <- 0
surveyuncon$republican[surveyuncon$pid7==3] <- -1
surveyuncon$republican[surveyuncon$pid7==2] <- -2
surveyuncon$republican[surveyuncon$pid7==1] <- -3


surveyuncon$faminc[surveyuncon$faminc==13] <- NA 
surveyuncon$faminc[surveyuncon$faminc==14] <- NA 


surveyuncon$biden2020 <- NA 
surveyuncon$biden2020[surveyuncon$voted_pres_20==1] <- 1 
surveyuncon$biden2020[surveyuncon$voted_pres_20==2] <- 0 


surveyuncon$presturnout20 <- NA 
surveyuncon$presturnout20[surveyuncon$voted_pres_20==1] <- 1 
surveyuncon$presturnout20[surveyuncon$voted_pres_20==2] <- 1 
surveyuncon$presturnout20[surveyuncon$voted_pres_20==3] <- 1 
surveyuncon$presturnout20[surveyuncon$voted_pres_20==4] <- 0 

surveyuncon$presturnout16 <- NA 
surveyuncon$presturnout16[surveyuncon$voted_pres_16==1] <- 1 
surveyuncon$presturnout16[surveyuncon$voted_pres_16==2] <- 1 
surveyuncon$presturnout16[surveyuncon$voted_pres_16==3] <- 1 
surveyuncon$presturnout16[surveyuncon$voted_pres_16==4] <- 0 


surveyuncon$presturnout <- NA 
surveyuncon$presturnout[surveyuncon$voted_pres_20==1] <- 1 
surveyuncon$presturnout[surveyuncon$voted_pres_20==2] <- 1 
surveyuncon$presturnout[surveyuncon$voted_pres_20==3] <- 1 
surveyuncon$presturnout[surveyuncon$voted_pres_20==4] <- 0 
surveyuncon$presturnout[surveyuncon$voted_pres_16==1] <- 1 
surveyuncon$presturnout[surveyuncon$voted_pres_16==2] <- 1 
surveyuncon$presturnout[surveyuncon$voted_pres_16==3] <- 1 
surveyuncon$presturnout[surveyuncon$voted_pres_16==4] <- 0 

surveyuncon$demsenatevote <- NA 
surveyuncon$demsenatevote[surveyuncon$voted_sen_party==1] <- 1
surveyuncon$demsenatevote[surveyuncon$voted_sen_party>=2] <- 0

surveyuncon$senatevote <- NA 
surveyuncon$senatevote[surveyuncon$voted_sen==1] <- 1
surveyuncon$senatevote[surveyuncon$voted_sen==2] <- 1
surveyuncon$senatevote[surveyuncon$voted_sen==3] <- 1
surveyuncon$senatevote[surveyuncon$voted_sen==4] <- 1
surveyuncon$senatevote[surveyuncon$voted_sen==5] <- 1
surveyuncon$senatevote[surveyuncon$voted_sen==6] <- 0
surveyuncon$senatevote[surveyuncon$voted_sen==8] <- 0

surveyuncon$demgovvote <- NA 
surveyuncon$demgovvote[surveyuncon$voted_gov_party==1] <- 1
surveyuncon$demgovvote[surveyuncon$voted_gov_party>=2] <- 0

surveyuncon$govvote <- NA 
surveyuncon$govvote[surveyuncon$voted_gov==1] <- 1
surveyuncon$govvote[surveyuncon$voted_gov==2] <- 1
surveyuncon$govvote[surveyuncon$voted_gov==3] <- 1
surveyuncon$govvote[surveyuncon$voted_gov==4] <- 1
surveyuncon$govvote[surveyuncon$voted_gov==5] <- 1
surveyuncon$govvote[surveyuncon$voted_gov==6] <- 0
surveyuncon$govvote[surveyuncon$voted_gov==8] <- 0

surveyuncon <- merge(surveyuncon, rural, by="FIPS")

# Proportion Covered Weight
surveyuncon$uncontesteddem2 <- (surveyuncon$uncontesteddem*surveyuncon$CountyPercent_updated)

# full county dummies for robustness
toy <- lm(contacted ~ uncontesteddem + republican + educ + faminc + female + black +  latinx + birthyr+ rural  + as.factor(Year) +as.factor(State.x), data=surveyuncon, w=weight_cumulative)

toy1 <- lm(biden2020 ~ uncontesteddem + republican + educ + faminc + female + black +  latinx + birthyr  + rural +  as.factor(Year) +as.factor(State.x), data=surveyuncon, w=weight_cumulative)

toy2 <- lm(demsenatevote ~ uncontesteddem + republican + educ + faminc + female + black +  latinx + birthyr + rural + as.factor(Year)+as.factor(State.x), data=surveyuncon, w=weight_cumulative)

toy3 <- lm(demgovvote ~ uncontesteddem + republican + educ + faminc + female + black +  latinx + birthyr + rural + as.factor(Year)+as.factor(State.x), data=surveyuncon, w=weight_cumulative)


toycl <- coeftest(toy, vcov = vcovCL, cluster = ~FIPS)
toy1cl <- coeftest(toy1, vcov = vcovCL, cluster = ~FIPS)
toy2cl <- coeftest(toy2, vcov = vcovCL, cluster = ~FIPS)
toy3cl <- coeftest(toy3, vcov = vcovCL, cluster = ~FIPS)

stargazer(toy, toy1, toy2, toy3, 
          covariate.labels=c("Uncontested Dem", "Republican", "Education", "Income", "Female", "Black", 
                             "Latinx", "Birth Year") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 

stargazer(toycl, toy1cl, toy2cl, toy3cl,
          covariate.labels=c("Uncontested Dem","Republican", "Education", "Income", "Female", "Black", 
                             "Latinx", "Birth Year") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 



# weighted by coverage (main text results)
toy <- lm(contacted ~ uncontesteddem2 + republican + educ + faminc + female + black +  latinx + birthyr+ rural  + as.factor(Year) +as.factor(State.x), data=surveyuncon, w=weight_cumulative)

toy1 <- lm(biden2020 ~ uncontesteddem2 + republican + educ + faminc + female + black +  latinx + birthyr  + rural +  as.factor(Year) +as.factor(State.x), data=surveyuncon, w=weight_cumulative)

toy2 <- lm(demsenatevote ~ uncontesteddem2 + republican + educ + faminc + female + black +  latinx + birthyr + rural + as.factor(Year)+as.factor(State.x), data=surveyuncon, w=weight_cumulative)

toy3 <- lm(demgovvote ~ uncontesteddem2 + republican + educ + faminc + female + black +  latinx + birthyr + rural + as.factor(Year)+as.factor(State.x), data=surveyuncon, w=weight_cumulative)



toycl <- coeftest(toy, vcov = vcovCL, cluster = ~FIPS)
toy1cl <- coeftest(toy1, vcov = vcovCL, cluster = ~FIPS)
toy2cl <- coeftest(toy2, vcov = vcovCL, cluster = ~FIPS)
toy3cl <- coeftest(toy3, vcov = vcovCL, cluster = ~FIPS)

stargazer(toy, toy1, toy2, toy3, 
          covariate.labels=c("Uncontested Dem", "Republican", "Education", "Income", "Female", "Black", 
                             "Latinx", "Birth Year") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 

stargazer(toycl, toy1cl, toy2cl, toy3cl,
          covariate.labels=c("Uncontested Dem","Republican", "Education", "Income", "Female", "Black", 
                             "Latinx", "Birth Year") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 


#############################################
# Texas Uncontested Seats County Results
############################################

# Loading Data
tx_cong <- read_csv("texas_congress_vacant.csv")

tx_shouse <- read_csv("texas_statehouse_vacant.csv")

tx_ssenate <- read_csv("texas_senate_vacant.csv")

ruraltx <- rural[which(rural$State=="tx"),]

# Merge Part I 
tx_cong <- merge(tx_cong, ruraltx, by.x = "County", by.y = "County_Name")
tx_shouse <- merge(tx_shouse, ruraltx, by.x = "County", by.y = "County_Name")
tx_ssenate <- merge(tx_ssenate, ruraltx, by.x = "County", by.y = "County_Name")

texas_electionsp2 <- texas_elections2

texas_electionsp2$Year <- texas_electionsp2$ELECTION_YEAR

# Merge Part II 
tx_congcounty <- merge(tx_cong, texas_electionsp2, by=c("FIPS", "Year"), all=T)
tx_shousecounty <- merge(tx_shouse, texas_electionsp2, by=c("FIPS", "Year"), all=T)
tx_ssenatecounty <- merge(tx_ssenate, texas_electionsp2, by=c("FIPS", "Year"), all=T)

tx_congcounty$area <- tolower(tx_congcounty$County.x)
tx_shousecounty$area <- tolower(tx_shousecounty$County.x)
tx_ssenatecounty$area <- tolower(tx_ssenatecounty$County.x)
txguber_2006_2022$area <- tolower(txguber_2006_2022$Area)

# Merge Part III 
tx_congcounty <- merge(tx_congcounty, txguber_2006_2022, by=c("area", "Year"), all=T)
tx_shousecounty <- merge(tx_shousecounty, txguber_2006_2022, by=c("area", "Year"), all=T)
tx_ssenatecounty <- merge(tx_ssenatecounty, txguber_2006_2022, by=c("area", "Year"), all=T)

# Coding Independent Variable
tx_congcounty$uncontesteddem_cong <- NA 
tx_congcounty$uncontesteddem_cong[tx_congcounty$`Democrat (Y/N)`=="Y"] <- 0 
tx_congcounty$uncontesteddem_cong[tx_congcounty$`Democrat (Y/N)`=="N"] <- 1

tx_congcounty$uncontestedgop_cong<- NA 
tx_congcounty$uncontestedgop_cong[tx_congcounty$`Republican (Y/N)`=="Y"] <- 0 
tx_congcounty$uncontestedgop_cong[tx_congcounty$`Republican (Y/N)`=="N"] <- 1

tx_shousecounty$uncontesteddem_sh <- NA 
tx_shousecounty$uncontesteddem_sh[tx_shousecounty$`Democrat (Y/N)`=="Y"] <- 0 
tx_shousecounty$uncontesteddem_sh[tx_shousecounty$`Democrat (Y/N)`=="N"] <- 1

tx_shousecounty$uncontestedgop_sh<- NA 
tx_shousecounty$uncontestedgop_sh[tx_shousecounty$`Republican (Y/N)`=="Y"] <- 0 
tx_shousecounty$uncontestedgop_sh[tx_shousecounty$`Republican (Y/N)`=="N"] <- 1

tx_ssenatecounty$uncontesteddem_ss <- NA 
tx_ssenatecounty$uncontesteddem_ss[tx_ssenatecounty$`Democrat (Y/N)`=="Y"] <- 0 
tx_ssenatecounty$uncontesteddem_ss[tx_ssenatecounty$`Democrat (Y/N)`=="N"] <- 1

tx_ssenatecounty$uncontestedgop_ss<- NA 
tx_ssenatecounty$uncontestedgop_ss[tx_ssenatecounty$`Republican (Y/N)`=="Y"] <- 0 
tx_ssenatecounty$uncontestedgop_ss[tx_ssenatecounty$`Republican (Y/N)`=="N"] <- 1

# Weighting Independent Variable by County Coverage
tx_congcounty$treatment_congdem <- (tx_congcounty$uncontesteddem_cong*tx_congcounty$CountyPercent)
tx_congcounty$treatment_conggop<- (tx_congcounty$uncontestedgop_cong*tx_congcounty$CountyPercent)
tx_shousecounty$treatment_shousedem <- (tx_shousecounty$uncontesteddem_sh*tx_shousecounty$CountyPercent)
tx_shousecounty$treatment_shousegop <- (tx_shousecounty$uncontestedgop_sh*tx_shousecounty$CountyPercent)
tx_ssenatecounty$treatment_ssenategop <- (tx_ssenatecounty$uncontestedgop_ss*tx_ssenatecounty$CountyPercent)
tx_ssenatecounty$treatment_ssenatedem <- (tx_ssenatecounty$uncontesteddem_ss*tx_ssenatecounty$CountyPercent)



keep1 <- c("FIPS", "Year", "area", "treatment_congdem")
keep2 <- c("FIPS", "Year", "area",  "treatment_shousedem")
keep3 <- c("FIPS", "Year", "area", "treatment_ssenatedem")

tx_congcounty2 <- tx_congcounty[keep1]
tx_shousecounty2 <- tx_shousecounty[keep2]
tx_ssenatecounty2  <- tx_ssenatecounty[keep3]


# Bringing all texas uncontested seats together post merge
tx_congcounty3 <- tx_congcounty2 %>%
  group_by(FIPS, Year) %>%
  slice(which.max(treatment_congdem))

tx_shousecounty3 <- tx_shousecounty2 %>%
  group_by(FIPS, Year) %>%
  slice(which.max(treatment_shousedem))

tx_ssenatecounty3 <- tx_ssenatecounty2 %>%
  group_by(FIPS, Year) %>%
  slice(which.max(treatment_ssenatedem))

texasuncontest <- merge(tx_congcounty3, tx_shousecounty3 , by=c("FIPS", "Year"), , all=T )
texasuncontest <- merge(texasuncontest, tx_ssenatecounty3 , by=c("FIPS", "Year"), , all=T )

keep <- c("FIPS", "Year", "area.x", "treatment_congdem", "treatment_shousedem", "treatment_ssenatedem")

texasuncontest <- texasuncontest[keep]

texasuncontest2 <- merge(texasuncontest, texas_electionsp2, by=c("FIPS", "Year"), all=T)

texasuncontest2$area <- texasuncontest2$area.x
# Remerging Elections Data
texasuncontest3 <- merge(texasuncontest2, txguber_2006_2022, by=c("area", "Year"), all=T)
# Remerging County Demo Data
texasuncontest3 <- merge(texasuncontest3, ruraltx, by="FIPS")

texasuncontest4 <- texasuncontest3[which(texasuncontest3$DemVotesMajorPercent>=0),]

texasuncontest4$CensusPop2 <- as.numeric(texasuncontest4$CensusPop)
texasuncontest4$govvote <- (texasuncontest4$DemVotesMajorPercent/100)

texasuncontest3$govvote <- (texasuncontest3$DemVotesMajorPercent/100)

texasuncontest3$Population_2010 <- texasuncontest3$Population_2010.x
texasuncontest4$Population_2010 <- texasuncontest4$Population_2010.y

# Presidential Vote 
toy1c <- lm(demshare ~ treatment_congdem + as.factor(Year), data=texasuncontest3, w=Population_2010)
toy2c <- lm(demshare ~ treatment_shousedem + as.factor(Year), data=texasuncontest3, w=Population_2010)
toy3c <- lm(demshare~ treatment_ssenatedem + as.factor(Year), data=texasuncontest3, w=Population_2010)
toy1 <- lm(demshare ~ treatment_congdem + as.factor(FIPS) + as.factor(Year), data=texasuncontest3, w=Population_2010)
toy2 <- lm(demshare ~ treatment_shousedem + as.factor(FIPS)+ as.factor(Year), data=texasuncontest3, w=Population_2010)
toy3 <- lm(demshare~ treatment_ssenatedem + as.factor(FIPS)+ as.factor(Year), data=texasuncontest3, w=Population_2010)


toy1cl <- coeftest(toy1, vcov = vcovCL, cluster = ~FIPS)
toy2cl <- coeftest(toy2, vcov = vcovCL, cluster = ~FIPS)
toy3cl <- coeftest(toy3, vcov = vcovCL, cluster = ~FIPS)

toy1clc <- coeftest(toy1c, vcov = vcovCL, cluster = ~FIPS)
toy2clc <- coeftest(toy2c, vcov = vcovCL, cluster = ~FIPS)
toy3clc <- coeftest(toy3c, vcov = vcovCL, cluster = ~FIPS)




stargazer(toy1c,toy1 ,toy2c, toy2, toy3c, toy3, 
          dep.var.labels=c("Democratic Presidential Vote"),
          covariate.labels=c("Congress", "State House", "State Senate") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 


stargazer(toy1clc, toy1cl, toy2clc, toy2cl, toy3clc, toy3cl,
          dep.var.labels=c("Democratic Presidential  Vote"),
          covariate.labels=c("Congress", "State House", "State Senate") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 


# Governor Vote
toy1c <- lm(govvote ~ treatment_congdem + as.factor(Year), data=texasuncontest4, w=Population_2010)
toy2c <- lm(govvote ~ treatment_shousedem + as.factor(Year), data=texasuncontest4, w=Population_2010)
toy3c <- lm(govvote~ treatment_ssenatedem + as.factor(Year), data=texasuncontest4, w=Population_2010)
toy1 <- lm(govvote ~ treatment_congdem + as.factor(FIPS) + as.factor(Year), data=texasuncontest4, w=Population_2010)
toy2 <- lm(govvote~ treatment_shousedem + as.factor(FIPS)+ as.factor(Year), data=texasuncontest4, w=Population_2010)
toy3 <- lm(govvote~ treatment_ssenatedem + as.factor(FIPS)+ as.factor(Year), data=texasuncontest4, w=Population_2010)


toy1cl <- coeftest(toy1, vcov = vcovCL, cluster = ~FIPS)
toy2cl <- coeftest(toy2, vcov = vcovCL, cluster = ~FIPS)
toy3cl <- coeftest(toy3, vcov = vcovCL, cluster = ~FIPS)

toy1clc <- coeftest(toy1c, vcov = vcovCL, cluster = ~FIPS)
toy2clc <- coeftest(toy2c, vcov = vcovCL, cluster = ~FIPS)
toy3clc <- coeftest(toy3c, vcov = vcovCL, cluster = ~FIPS)



stargazer(toy1,toy2,   toy3, 
          dep.var.labels=c("Democratic Gubernatorial Vote"),
          covariate.labels=c("Congress", "State House", "State Senate") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 


stargazer(toy1cl, toy2cl, toy3cl,
          dep.var.labels=c("Democratic Gubernatorial Vote"),
          covariate.labels=c("Congress", "State House", "State Senate") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 



#############################################
# Texas Uncontested Seats CES Results
############################################

CESTEXAS <- read.csv("cestexas.csv")
CESTEXAS2 <- CESTEXAS
CESTEXAS2$FIPS <- CESTEXAS2$FIPS.x
txuncontextces <- merge(texasuncontest3, CESTEXAS2, by=c("FIPS", "Year"))

# Coding DVs
txuncontextces$demvote_gov <- NA 
txuncontextces$demvote_gov[txuncontextces$voted_gov==1] <- 1 
txuncontextces$demvote_gov[txuncontextces$voted_gov==2] <- 0 

txuncontextces$dempresvote <- NA 
txuncontextces$dempresvote[txuncontextces$voted_pres_party==1] <- 1 
txuncontextces$dempresvote[txuncontextces$voted_pres_party==2] <- 0
txuncontextces$dempresvote[txuncontextces$voted_pres_party==3] <- 0

txuncontextces$votedemsenate <- NA 
txuncontextces$votedemsenate[txuncontextces$voted_sen_party==1] <- 1
txuncontextces$votedemsenate[txuncontextces$voted_sen_party>=2] <- 0

# Coding Demographics
txuncontextces$black <- 0
txuncontextces$black[txuncontextces$race==2] <- 1

txuncontextces$latino <- 0
txuncontextces$latino[txuncontextces$race==3] <- 1

txuncontextces$faminc[txuncontextces$faminc==13] <- NA

txuncontextces$republican <- NA 
txuncontextces$republican[txuncontextces$pid7==1] <- -3
txuncontextces$republican[txuncontextces$pid7==2] <- -2
txuncontextces$republican[txuncontextces$pid7==3] <- -1
txuncontextces$republican[txuncontextces$pid7==4] <- 0
txuncontextces$republican[txuncontextces$pid7==5] <-1
txuncontextces$republican[txuncontextces$pid7==6] <-2
txuncontextces$republican[txuncontextces$pid7==7] <-3

txuncontextces$female <- 0
txuncontextces$female[txuncontextces$gender==2] <-1 


# Pres vote 
toy3c <- lm(dempresvote ~ treatment_congdem  + republican + educ + faminc + black + latino + birthyr + as.factor(Year), data= txuncontextces, w=weight_cumulative)
toy3cl <- coeftest(toy3c, vcov = vcovCL, cluster = ~FIPS)
toy3s <- lm(dempresvote ~ treatment_ssenatedem  + republican + educ + faminc + black + latino + birthyr + as.factor(Year), data= txuncontextces, w=weight_cumulative)
toy3scl <- coeftest(toy3s, vcov = vcovCL, cluster = ~FIPS)
toy3h <- lm(dempresvote ~ treatment_shousedem + republican + educ + faminc + black + latino + birthyr + as.factor(Year), data= txuncontextces, w=weight_cumulative)
toy3hcl <- coeftest(toy3h, vcov = vcovCL, cluster = ~FIPS)


stargazer(toy3c, toy3s, toy3h,
          dep.var.labels=c("Democratic Presidential Vote"),
          covariate.labels=c("Congress", "State House", "State Senate") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 


stargazer(toy3cl, toy3scl, toy3hcl,
          dep.var.labels=c("Democratic Gubernatorial Vote"),
          covariate.labels=c("Congress", "State House", "State Senate") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 



toy3c <- lm(dempresvote ~ treatment_congdem  + republican + educ + faminc + black + latino + birthyr + as.factor(Year) +as.factor(FIPS), data= txuncontextces, w=weight_cumulative)
toy3cl <- coeftest(toy3c, vcov = vcovCL, cluster = ~FIPS)
toy3s <- lm(dempresvote ~ treatment_ssenatedem  + republican + educ + faminc + black + latino + birthyr + as.factor(Year)+as.factor(FIPS), data= txuncontextces, w=weight_cumulative)
toy3scl <- coeftest(toy3s, vcov = vcovCL, cluster = ~FIPS)
toy3h <- lm(dempresvote ~ treatment_shousedem + republican + educ + faminc + black + latino + birthyr + as.factor(Year)+as.factor(FIPS), data= txuncontextces, w=weight_cumulative)
toy3hcl <- coeftest(toy3h, vcov = vcovCL, cluster = ~FIPS)

stargazer(toy3c, toy3s, toy3h,
          dep.var.labels=c("Democratic Presidential Vote"),
          covariate.labels=c("Congress", "State House", "State Senate") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 


stargazer(toy3cl, toy3scl, toy3hcl,
          dep.var.labels=c("Democratic Gubernatorial Vote"),
          covariate.labels=c("Congress", "State House", "State Senate") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 


# Gov vote 
toy4c <- lm(demvote_gov ~ treatment_congdem  + republican + educ + faminc + black + latino + birthyr + as.factor(Year), data= txuncontextces, w=weight_cumulative)
toy4cl <- coeftest(toy4c, vcov = vcovCL, cluster = ~FIPS)
toy4s <- lm(demvote_gov ~ treatment_ssenatedem  + republican + educ + faminc + black + latino + birthyr + as.factor(Year), data= txuncontextces, w=weight_cumulative)
toy4scl <- coeftest(toy4s, vcov = vcovCL, cluster = ~FIPS)
toy4h <- lm(demvote_gov ~ treatment_shousedem + republican + educ + faminc + black + latino + birthyr + as.factor(Year), data= txuncontextces, w=weight_cumulative)
toy4hcl <- coeftest(toy4h, vcov = vcovCL, cluster = ~FIPS)

stargazer(toy4c, toy4s, toy4h,
          dep.var.labels=c("Democratic Presidential Vote"),
          covariate.labels=c("Congress", "State House", "State Senate") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 


stargazer(toy4cl, toy4scl, toy4hcl,
          dep.var.labels=c("Democratic Gubernatorial Vote"),
          covariate.labels=c("Congress", "State House", "State Senate") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 

toy4c <- lm(demvote_gov ~ treatment_congdem  + republican + educ + faminc + black + latino + birthyr + as.factor(Year)+as.factor(FIPS), data= txuncontextces, w=weight_cumulative)
toy4cl <- coeftest(toy4c, vcov = vcovCL, cluster = ~FIPS)
toy4s <- lm(demvote_gov ~ treatment_ssenatedem  + republican + educ + faminc + black + latino + birthyr + as.factor(Year)+as.factor(FIPS), data= txuncontextces, w=weight_cumulative)
toy4scl <- coeftest(toy4s, vcov = vcovCL, cluster = ~FIPS)
toy4h <- lm(demvote_gov ~ treatment_shousedem + republican + educ + faminc + black + latino + birthyr + as.factor(Year)+as.factor(FIPS), data= txuncontextces, w=weight_cumulative)
toy4hcl <- coeftest(toy4h, vcov = vcovCL, cluster = ~FIPS)

stargazer(toy4c, toy4s, toy4h,
          dep.var.labels=c("Democratic Presidential Vote"),
          covariate.labels=c("Congress", "State House", "State Senate") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 


stargazer(toy4cl, toy4scl, toy4hcl,
          dep.var.labels=c("Democratic Gubernatorial Vote"),
          covariate.labels=c("Congress", "State House", "State Senate") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 

# Senate  Dem 
toy6c <- lm(votedemsenate~ treatment_congdem  + republican + educ + faminc + black + latino + birthyr + as.factor(Year), data= txuncontextces, w=weight_cumulative)
toy6cl <- coeftest(toy6c, vcov = vcovCL, cluster = ~FIPS)
toy6s <- lm(votedemsenate ~ treatment_ssenatedem  + republican + educ + faminc + black + latino + birthyr + as.factor(Year), data= txuncontextces, w=weight_cumulative)
toy6scl <- coeftest(toy6s, vcov = vcovCL, cluster = ~FIPS)
toy6h <- lm(votedemsenate ~ treatment_shousedem + republican + educ + faminc + black + latino + birthyr + as.factor(Year), data= txuncontextces, w=weight_cumulative)
toy6hcl <- coeftest(toy6h, vcov = vcovCL, cluster = ~FIPS)

stargazer(toy6c, toy6s, toy6h,
          dep.var.labels=c("Democratic Presidential Vote"),
          covariate.labels=c("Congress", "State House", "State Senate") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 


stargazer(toy6cl, toy6scl, toy6hcl,
          dep.var.labels=c("Democratic Gubernatorial Vote"),
          covariate.labels=c("Congress", "State House", "State Senate") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 

toy6c <- lm(votedemsenate~ treatment_congdem  + republican + educ + faminc + black + latino + birthyr + as.factor(Year)+as.factor(FIPS), data= txuncontextces, w=weight_cumulative)
toy6cl <- coeftest(toy6c, vcov = vcovCL, cluster = ~FIPS)
toy6s <- lm(votedemsenate ~ treatment_ssenatedem  + republican + educ + faminc + black + latino + birthyr + as.factor(Year)+as.factor(FIPS), data= txuncontextces, w=weight_cumulative)
toy6scl <- coeftest(toy6s, vcov = vcovCL, cluster = ~FIPS)
toy6h <- lm(votedemsenate ~ treatment_shousedem + republican + educ + faminc + black + latino + birthyr + as.factor(Year)+as.factor(FIPS), data= txuncontextces, w=weight_cumulative)
toy6hcl <- coeftest(toy6h, vcov = vcovCL, cluster = ~FIPS)


stargazer(toy6c, toy6s, toy6h,
          dep.var.labels=c("Democratic Presidential Vote"),
          covariate.labels=c("Congress", "State House", "State Senate") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 


stargazer(toy6cl, toy6scl, toy6hcl,
          dep.var.labels=c("Democratic Gubernatorial Vote"),
          covariate.labels=c("Congress", "State House", "State Senate") , style = "ajps", 
          omit = c("factor"),
          omit.stat=c("LL","ser","f"), no.space=FALSE, align=TRUE) 





#############################################
# # SAVED OUTPUT FOR COEFFICIENT PLOTS 
############################################
texaseffect  <- read_csv("texas_uncontest_effects.csv")


limits = aes(ymax = Effect + (1.64*SE), ymin=Effect- (1.64*SE))

dodge = position_dodge(width=0.9)

apatheme=theme_bw()+
  theme(panel.grid.major=element_blank(),
        panel.grid.minor=element_blank(),
        panel.border=element_blank(),
        axis.line=element_line(),
        text=element_text(family='Times'))

demsen <- texaseffect[which(texaseffect$Outcome=="Senate"),]
dempres <- texaseffect[which(texaseffect$Outcome=="President"),]
demgov<- texaseffect[which(texaseffect$Outcome=="Governor"),]



p=ggplot(dempres, aes(x = Variable, y = Effect, color=Model, fill=Model))+
  geom_bar(stat='identity', position=dodge)+
  geom_errorbar(limits, position=dodge, width=0.25, , color="black")+
  apatheme + geom_hline(yintercept = 0, linetype="dashed", color="black") + 
  labs(y="Dem Presidential Vote", color="Model", fill="Model", x="Uncontested Seat Type") + theme(text = element_text(size = 17)) + 
  theme(legend.position="bottom") + scale_color_grey() + scale_fill_grey()
p

p2=ggplot(demsen, aes(x = Variable, y = Effect, color=Model, fill=Model))+
  geom_bar(stat='identity', position=dodge)+
  geom_errorbar(limits, position=dodge, width=0.25, , color="black")+
  apatheme + geom_hline(yintercept = 0, linetype="dashed", color="black") + 
  labs(y="Dem Senate Vote", color="Model", fill="Model", x="Uncontested Seat Type")  + theme(text = element_text(size = 17)) + 
  theme(legend.position="bottom") + scale_color_grey() + scale_fill_grey()
p2


p3=ggplot(demgov, aes(x = Variable, y = Effect, color=Model, fill=Model))+
  geom_bar(stat='identity', position=dodge)+
  geom_errorbar(limits, position=dodge, width=0.25, , color="black")+
  apatheme + geom_hline(yintercept = 0, linetype="dashed", color="black") + 
  labs(y="Dem Governor Vote", color="Model", fill="Model", x="Uncontested Seat Type")  + theme(text = element_text(size = 17)) + 
  theme(legend.position="bottom") + scale_color_grey() + scale_fill_grey()
p3

library(ggpubr)

blank <- ggplot() + theme_void()


ggarrange(
  ggarrange(p, p2, blank, nrow = 1,
            widths = c(1, 1, 0), legend = "none"),
  ggarrange(blank, p3, blank, nrow = 1,
            widths = c(0.5, 1, 0.5), legend = "bottom"),
  nrow = 2,
  legend = "bottom", common.legend = T
)

ggarrange(p, p2, p3, nrow=2, ncol=2, common.legend = T, legend = "bottom")


